We are IntechOpen, 
the world’s leading publisher of 


Open Access books 
Built by scientists, for scientists 


5.300 130,000 150M 


ailable International authors and editors Downloads 


Our author among the 


154 TOP 1% 12.2% 


Countries delivered to most cited s Contributors from top 500 universities 


Selection of our books indexed in the Book Citation Index 
in Web of Science™ Core Collection (BKCI) 


Interested in publishing with us? 
Contact book.department@intechopen.com 


Numbers displayed above are based on latest data collected. 
F information visit www.intechopen.com 


ay 


Chapter 


Topological Visualisation 
Techniques for Volume Multifield 
Data 


Dean P. Thomas, Rita Borgo, Robert S. Laramee 
and Simon J. Hands 


Abstract 


This survey paper provides an overview of topological visualisation techniques 
for scalar data sets. Topological algorithms are used to reduce scalar fields to a 
skeleton by mapping critical changes in the topology to the vertices of graph struc- 
tures. These can be visualised using graph drawing techniques or used as a method 
of seeding meshes of distinct objects existing in the data. Many techniques are 
discussed in detail, beginning with a review of algorithms working on scalar fields 
defined with a single variable, and then generalised to multivariate and temporal 
data. The survey is completed with a discussion of methods of presenting data in 
higher dimensions. 


Keywords: indirect volume rendering, topology driven visualisation, 
multivariate visualisation, contour tree, Reeb graph, Reeb space, joint contour net, 
Reeb skeleton, scalar data 


1. Introduction 


The Marching Cubes (MC) algorithm [1] is a long established method used in 
indirect visualisation for creating mathematical models of data existing in a scalar 
field. Early research centred upon improvements for the presentation of data from 
medical sources such as CT and MRI scans. Whilst well suited for approximating 
isosurfaces on smooth functions, MC derived algorithms typically struggle to accu- 
rately capture features such as sharp edges and corners. These features can often to 
lead to poorly shaped triangles. Various approaches have been suggested for correct- 
ing these limitations including Extended Marching Cubes [2], Dual Marching Cubes 
[3], and Cube Merging [4]. The algorithm continues to be refined and improved for 
various purposes, and remains an active field of study [5]. 

Besides improvements to model quality, a recurring theme in research linked 
to the algorithm is identifying and solving topological ambiguities in the models 
it creates. Typical ambiguities result in the absence of triangles between two 
surfaces existing in opposing corners of a MC cell. This can lead to creation of a 
model which contains two separate objects, when in reality the model should be a 
single joined object. Many approaches have been suggested for overcoming these 
ambiguities including the use of Marching Tetrahedra [6]. The use of topological 
algorithms removes the possibility for ambiguities, whilst providing a method for 
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seeding disjoint contours using similar lookup tables to those of Marching Cubes. 
Algorithms for surface generation, such as Marching Cubes, can also be generalised 
for use in multivariate data sets. 

The remainder of this chapter is structured as follows. Section 2 describes topo- 
logical visualisation techniques for scalar fields that solve the ambiguity problems 
of MC. In Section 3 we discuss extending topological techniques to cover multivari- 
ate data sets. Section 4 considers displaying the output of topological algorithms 
from higher dimensional data sets. The chapter is concluded in Section 5. 


2. Topology driven visualisation 


Algorithms for computing isosurfaces, otherwise known as level sets, are 
prone to topological inaccuracies [7-9]. During the construction of an isosurface 
no topological information is used; hence, an isosurface is unable to distinguish 
individual connected regions. Morse theory provides a method for computing 
the topology of scalar data, as we sweep through the isovalue range, by consid- 
ering the gradient and second derivatives of the level set at each vertex in the 
scalar field. Through the use of tetrahedral mesh cells (for 3D scalar fields) we 
can approximate the isosurface as a piecewise linear function between sampling 
points. However, at the boundary between cells this approach invalidates a key 
condition of Morse theory which requires a function to have derivatives defined 
at all sampling points. Work by Edelsbrunner et al. [10] and Bremer et al. [11] 
provides a mechanism for handling these limitations, allowing us to consider the 
data as a Morse function. 

Topology driven visualisation is used to capture characteristics of a data set 
using basic topological invariants such as the number of holes or points in a 
connected region. This information provides a means for obtaining a skeleton of 
the data which can then be used to construct visualisations in various forms. Use 
of topological visualisation techniques can see a reduction in the size of a large 
data set, as only key features need to be retained. By using topological structures 
to represent a data set intuitive methods can be used for speeding up computa- 
tions; for example, by allowing many parts of a data set to be rendered in parallel. 
A further increase in efficiency is provided by the ability to bypass redundant 
computations—such as the processing of empty cells in MC. Increases in rendering 
performance come with the overhead of the need to do a one-time pre-processing 
of the data. 

Many topology based algorithms rely on abstract graph based structures for 
storage of data. In many situations the graph structures can be viewed directly to 
form an overview of the data in an easily perceived format. For those with little 
prior knowledge of topological visualisation techniques it can be hard to form an 
initial understanding of the data using such formats. In order to better understand 
the data it is common to directly link the graph visualisation, either statically or 
dynamically, with a rendered view of the data. However, difficulty in understand- 
ing topological data representations can be further complicated when a data set is 
large in size. Hence, other approaches are sometimes used that try to take advantage 
of human perception, such as topological landscapes. 


2.1 Terminology 


In order to present the algorithms and structures used in topological visualisa- 
tion the following key terminology is used. 
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2.1.1 Simplex 


The simplex is a generalisation of the triangle or tetrahedron to n dimensions. 
Properties of simplexes in dimensions 0 (a point) through to 4 (a Pentatope) are 
given in Table 1. Computational topology often refers to meshes as simplicial 
complexes—a set of n simplexes glued together at their boundaries. 


2.1.2 Betti numbers 


Enable us to categorise the topology of a mesh existing in n dimensions by exam- 
ining the connectivity of its structure as a number of n dimensional simplexes. Each 
dimension, from o ton - 1, has its own Betti number f,, associated with it to represent 
the quantity of cuts required through a mesh in order to separate it into two parts. 

We can consider these numbers as representing the number of holes in the mesh in n 
dimensions (Table 2). Most topological visualisation algorithms consider only the .oth 
dimensional Betti number, used to represent connectivity. However, it is possible to 
augment the graphs with higher order Betti numbers to provide feedback on changes 
in the mesh such as the morphing from a sphere to a torus (Table 3). In this work we 
are primarily only interested in the zeroth dimensional Betti number $o. 


2.1.3 Genus 


This is a way of expressing the number of holes in a surface. In three dimensions, 
a sphere is genus 0, and torus genus 1. Surfaces with additional holes are commonly 
referred to as an n-torus. 


2.1.4 Critical point 


The sampling points in a scalar field fœ), at which the function value is defined, 
can be categorised into two types. Critical points have a local gradient fx) =- o and 
represent positions where the topology of the level set changes. These relate to local 
extrema or saddle points in the data; the second derivative f (x) allows us to classify 


Name Dimension Common name Vertices Edges Faces Cells 
0-simplex 0 Point 1 — — — 
1-simplex 1 Line-segment 2 1 — — 
2-simplex 2 Triangle 3 3 1 — 
3-simplex 3 Tetrahedron 4 6 4 1 
4-simplex 4 Pentatope (or 5-cell) 5 10 10 5 
Table 1. 
Classification of simplexes in dimensions zero to four. 
Betti number Associated simplex Property 
Bo Point Connected components 
Py Line-segment Holes 
Po Triangles Voids 
Table 2. 


Betti numbers ba Pa and p, and associated concepts of connectivity. 
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Topological concept Bo By b2 

Sphere 1 0 1 

Torus 1 2 1 
Table 3. 


Betti numbers used to describe the sphere and torus. 


if the critical point is a local minimaf’(«) > 0 or a local maximaf’(x) < o. In the case 
where a critical point is a local minima all neighbouring vertices will have a higher 
function value, whilst a local maxima means that all surrounding vertices are lower 
in value. Regular points have a non-zero gradient and have no overall effect on the 
connectivity the level set. 


2.1.5 Contour 


The boundary of connected regions of a level set, also known as the 0-dimensional 
homology group, is known as contours. Each individual contour is an element in the 
level set that represents a distinct topological object within the isosurface. 


2.1.6 Topological persistence 


The most simplistic way of defining topological persistence is as a quantitative 
measure of importance of each contour. In literature topological persistence is 
often known as geometric measures [12] due to the fact it relates to measures such 
as volume or surface area. It is possible to compute persistence values directly from 
the meshes representing each contour by performing a piecewise sum of all cells in 
a connected region. Alternatively, the persistence can be approximated by counting 
the number of regular vertices making up a connected region. Persistence measures 
are often used as a method of eliminating noise from the topological structure by 
iteratively removing the contour with smallest persistence value. 

Most of the techniques and structures discussed in this section can be gener- 
alised for use in any number of dimensions—a list of general terms is provided in 
Table 4. 


2.2 Merge trees 


One of the fundamental data structures in computational topology are merge 
trees, commonly used to describe two similar concepts. Merge trees track connected 
level sets within the data as the isovalue is swept across the scalar range. The join tree 
captures connected sub-level sets, or regions, of the data below a specific isovalue 
threshold. Similarly the split tree captures connected super-level sets. 


Name Dimension Common name 

0-manifold 0 Point 

1-manifold 1 Polyline 

2-manifold 2 Polygon 

3-manifold 3 Polyhedron 

4-manifold 4 4-Polytope 
Table 4. 


Classification of manifolds in dimensions zero to four. 
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Merge trees present a simple method for capturing topological information as 
the arcs in a tree-like structure. This can then be used to compute zero-dimensional 
persistence measures representing the number of connected components in a region 
[14]. Merge trees are often used as a computation step in the construction of more 
complex topological data structures such as the contour tree, described in detail in 
Section 2.3. Alternatively, merge trees can be used directly for analysis and feature 
extraction in complex data sets [13]. 

A visual metaphor for understanding the concept of merge trees is given in 
Figure 1. To capture the merge tree, in this case the join tree of the grey func- 
tion (top-left), the isovalue is gradually decreased from the global maxima. In 
the top-centre image this has revealed three local maxima, each denoted as a red 
critical point. In the top right image a fourth local maxima has been revealed (the 
green contour) and the red and blue structures have merged, this is marked by 
a node in the graph. In the bottom left image, the pink and green surfaces have 
merged to form the cyan region. The bottom centre image sees the small intervals 
represented by the cyan and yellow surfaces merged into one large purple region. 
If the isovalue continues to be reduced it is possible to observe no further changes 
in the topology until reaching the global minima represented by the blue graph 
node. 


2.3 The contour tree 


The contour tree, introduced by Boyell and Ruston [16] uses concepts from 
Morse theory [17, 18] to capture changes in the topology of scalar field at critical 
points in the data. The tree structure allows the tracking of splits and joins in the 
topology as the isovalue is varied by tracing a path through the graph. Figure 2 
demonstrates how individual contours in the data are captured by the contour 
tree using a one-to-one correspondence with edges. Each leaf in the contour tree 
represents an individual local extrema and critical vertices represent a change in the 
connectivity of the level set. Typically, the contour tree does not contain informa- 
tion regarding changes in topological genus; however, this information can be 
added to the output [19]. 

We can formally describe the structure of the contour tree as follows: 


oli 


Figure 1. 

A visual representation of the merge tree algorithm for the function shown in the top-left. The merge tree 
captures each of the local maxima as the function height is reduced (moving left-to-right, top-to-bottom). A 
real world analogy is to think of the isovalue as the water level in a valley that gradually drains to reveal the 
peaks of multiple hills. Image courtesy of Landge et al. [13]. 
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Figure 2. 
Left: a surface that varies as the function value is varied. Right: a labelled contour tree representation of the 
surface. Image courtesy of Carr et al. [15]. 


e Vertices or supernodes 


o Leaf nodes correspond to: 
a. Local maxima where a contour is created. 


b.Local minima where a contour is destroyed. 


o Interior nodes correspond to: 


a. A saddle point where one more contours are created as a contour splits 
into two or more disjoint contours. 


b.A saddle point where one or more contours are destroyed as two or more 
disjoint contours merge into one. 


e Edges or superarcs 


o Represent a single contour between the supernode where it is created and 
the supernode where it is destroyed. 


In addition to optimizations of contour tree computation, research is focused 
on approaches for storing the contour tree hierarchy, allowing it to be traversed in 
an efficient manner. A number of related approaches have been proposed including 
Kd-trees, segment trees, and interval trees. These formats prove problematic as they 
can require a large amount of storage and are difficult to navigate. An alternative 
method is proposed in [20] using the observation that an isosurface is the level set of 
a continuous function in 3D space; therefore, a whole contour can be traced froma 
single element. Vertices where contours can be built from are given the name seeds; 
algorithms for computing seeds relate to work in image processing requiring similar 
storage facilities. An entire contour can be computed for the specified desired isovalue 
from the seed that is bounded by the two supernodes representing the critical points 
of the contour. Results showed that by using this method the number of seeds required 
for representation were in most cases significantly reduced. However, often the overall 
storage size of the contour tree was increased in comparison to other methods. 

Carr et al. [15] investigated the use of contour trees in higher dimensional 
data sets, whilst also improving upon the algorithm proposed in [21]. This also 
introduced the concept of augmented contour trees, an extension that added 
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non-branching vertices at non-critical points in the data to provide additional 
values for isosurfaces to be seeded from. This feature was built into a GUI allow- 
ing the user to identify regions of interest, using colour coding and to distinguish 
them as directed by the user [22]. Carr and Snoeyink [23] use the contour tree as 
the underlying data structure to generate object meshes using path seeds. The use 
of the contour tree can be extended to finding paths between two points given 
condition clauses, such as a minimum or maximum values, on a function defining a 
landscape [24]. 

Chiang et al. [25] introduce a modified contour tree construction algorithm 
that improves processing time by only sorting critical vertices in the merge tree 
construction stage. This vastly increases processing speeds in very large data sets. 
However, as observed by the authors, critical vertices are difficult to identify in 
data with dimensionality of four or more. Further increases in speed are offered 
by storing multiple seeds for the monotone path construction algorithm [26] 
used to generate surfaces. An increased storage overhead is required; however, a 
surface does not require complete re-extraction each time the isovalue is changed. 
Recently Gueunet et al. [27] presented a multithreaded approach to computing the 
augmented contour tree using a shared memory approach. Initial evaluations of 
the “contour forest” algorithm showed quicker computations achievable using the 
approach; however, at present the extent of the speed up is limited by load imbal- 
ance and redundant computations. 

For a given technique of subdividing the input domain, the output of the con- 
tour tree algorithm is fixed. However, in the case that a different technique is used to 
define the neighbourhood of sampling points the location of the critical point may 
change trivially. For example, changing the method of defining neighbours during 
split/merge tree computations can slightly alter the location of critical points. The 
underlying structure of the tree, the number of supernodes and superarcs, remains 
fixed for a given input with only the geometric locations of supernodes changing. 


2.4 The Reeb graph 


Using the contour tree comes with the limitation that the input must be defined 
on asimply connected domain that is free of loops. Limitations in the merge tree com- 
putation phase of the contour tree algorithm prevent it from correctly handling data 
without a boundary. However, the Reeb graph [28], a generalisation of the contour 
tree, can be used to compute the topology of scalar data in such situations [29]. As 
with the contour tree, the Reeb graph can be applied to models of any dimension 
provided it is represented on a simplicial mesh [30]. The Reeb graph has also been 
used to assist in the design of transfer functions in volume visualisation [31], by 
assigning opacity based upon how many objects were obscured by nested surfaces 
and their proximity to the edge of a scalar field. 

The first use of the Reeb graph for encoding topological features for visualisation 
purposes was by Shinagawa and Kunii [32] where the structure was used as a way of 
representing objects obtained from computerised-tomography (CT) sources [33]. 
An online algorithm is given by Pascucci et al. [30] that allowed the Reeb graph of 
a function to be computed using a streaming approach with the output continually 
updated as additional data points are added. Efficiency of the algorithm is optimal 
for two dimensional inputs and the streaming nature of the algorithm limits peak 
memory usage. However, the oln’) complexity of the algorithm, where n is the num- 
ber of triangles, means it performs less favourably for higher dimensional inputs. 

A key feature of many Reeb graph algorithms is that a 2-skeleton of the input 
volume is used to define the relationship between vertices, edges and triangles. 

The same restriction applies to scalar fields in n dimensions; hence, provided a 
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triangulation of the data is available the algorithm will be able to compute a cor- 
rect n dimensional Reeb graph. Reeb graph algorithms can largely be split into two 
groups; those that are sweep based, requiring the maintenance of level sets, and 
those that use a split and compute method to collapse the problem to that of a contour 
tree computation. A third alternative was given by Harvey et al. [34] that randomly 
collapsed triangles for arbitrary 2-skeletons to improve the running time to O(m logm). 

The sensitivity of Reeb graph computations to the topological genus of the input 
domain was first demonstrated by McLaughlin [29] in the case of a non-simply 
connected domain, such as that representing a sphere or torus. Additionally, the 
number of loops present in the data as a result of the Morse function can further 
extend the complexity of the algorithm. In order to address this undesired effect 
the Reeb graph computation can be reduced to that of the simpler contour tree 
algorithm. This approach was first used by Tierny et al. [35] where they performed 
“loop-surgery” by making symbolic cuts during the merge tree step of the algorithm. 
The symbolic cuts could then be stitched back together to retrieve a topologically 
correct Reeb graph of the data. An alternative approach, capable of generalising 
to higher dimensional inputs, was later proposed by Doraiswamy and Natarajan 
[36] that explicitly maintained level sets, eliminating the need for a loop-surgery 
pre-processing step. 

More recently Doraiswamy and Natarajan [37] offered an improved algorithm 
for constructing the Reeb graph of n-dimensional scalar fields, reverting to the 
split and compute contour tree technique, as originally proposed in [35]. A further 
optimised algorithm [38], using a variation of the split and compute technique, 
used the join tree of the data to identify potential loops. An important modification 
upon the technique developed by Tierny et al. [35] was the segmentation of domain 
into multiple loop free contour trees, instead of a single contour tree, thus enabling 
multiple regions of the input to be computed in parallel. 

Output from Reeb graph algorithms remains fixed provided the simplicial 
subdivision defined upon the input domain does not change. As is the case with the 
contour tree, variations can slightly alter the location of the critical vertices without 
changing the overall structure of the graph. 


2.5 Seeded contours 


Contours, as opposed to isosurfaces, can be grown using values extracted from 
the scalar field to produce topologically correct meshes. De Berg and van Kreveld 
[24] made the observation that a whole contour could be traced starting from a 
single seed element. Wyvill et al. [39] generate contour surfaces as a two-stage 
process; first, the core region is flood filled by evaluating neighbouring cells for 
inclusion in the level set. In the second stage boundaries are calculated using look- 
up tables similar in form to those of the MC algorithm [1]. 

Contour trees are a reliable source of seed cells for propagation algorithms similar 
to those of Wyvill et al. [39]; hence, the contour tree can be used to generate surfaces 
for each superarc at a given isovalue. Carr and Snoeyink [23] use the contour tree 
as the underlying data structure to store and generate object meshes. Rather than 
using minimal seed sets as used by van Kreveld et al. [20] to generate contours, an 
alternative method was deployed using path seeds. Carr et al. [15] also investigated 
the use of contour trees in higher dimensional datasets, whilst also improving upon 
the algorithm proposed by Tarasov and Vyalyi [21]. This enables users to identify 
individual contours of interest and distinguish them accordingly [22]. 

The flexible isosurface [40] is a technique for combining many of the features 
discussed in this section with the concepts of topological persistence (see Section 2.6) 
into a single interface. Simplification methods, as discussed in [12], can also be 
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applied to the data so as to display isosurfaces in a meaningful way. This approach 
allows each contour to be hidden or displayed according to the user preferences at 
runtime. In order to further aid data exploration contours can also remain fixed 
in view as the global isovalue is varied. Colour can be applied to each surface (ora 
sub-tree of the main contour tree) to allow assignment of meaning to contours by 
the user by providing a simple grouping mechanism. 


2.6 Topological persistence and simplification 


Different scientific fields of study often define the interesting features and 
attributes in a domain specific form. However, some features are invariant and 
can be useful in any field of science. The contour spectrum [41] was introduced 
as a method of relaying quantitative information about individual contours in 
scalar data including surface area and volume. Carr et al. [42] directly compare 
isosurface statistics against raw histograms of scalar data for a number of data sets. 
Measurements evaluated included the cell intersection count, triangle count and 
isosurface area. It was found that using these measures a truer distribution of the 
scalar field could be computed. An improvement was given by Meyer et al. [43] 
using concepts from geometric measure theory that minimised the effect of noise 
on the observed distributions. The key to this improvement was introducing a 
normalisation of the individual contour statistics to the domain average. 

Scalar field data often contains noise; when partitioned using a topology sensi- 
tive algorithm this manifests in the generation of a large quantity of small objects 
present at a limited range of isovalues. Methods for reducing noise were first 
proposed by Edelsbrunner et al. [14] using an iterative process that performed 
topological simplification by assessing the Betti-numbers of topological objects. The 
goal of the algorithm is to simplify shape whilst preserving the underlying topologi- 
cal features of data existing on a triangular mesh. 

Carr et al. [12] proposed the use of concepts originally discussed as part of the 
contour spectrum [41], such as enclosed volume and surface area, as an aid for noise 
removal. Objects are queued in ascending order, according to the user selected mea- 
sure, and iteratively removed until a terminating level of simplification is achieved. 
The effects of three different simplification features are compared using X-ray data 
from a human skull, where it was found that use of the isovalue range persistence 
measure (Figure 3) could result in the removal of important features such as the 
eyes. Carr et al. remark that simplification techniques are sensitive to the source 
and should be chosen using domain specific knowledge. The technique is suited to n 
dimensional data; however, time variate data requires additional considerations in 
the simplification process due to connectivity between time steps. 


2.7 Topology in direct volume rendering 


The techniques discussed in Section 2.5 model topological objects as meshes, 
meaning they are well suited to indirect volume visualisation. However, topology 
based approaches can also be applied to data displayed using direct volume render- 
ing approaches. 

A segmentation algorithm, such as the contour tree or volume skeleton tree, is first 
used to look for boundaries between objects in the scalar field topology. Following 
this, traditional direct volume rendering techniques can be applied to the data based 
upon attributes of the identified objects. Takeshima et al. [44] use attributes such 
as the number of equal valued contours and occlusion to assign levels of opacity to 
objects. The net effect was that the outermost objects were assigned lower opacities 
so as to not obscure features centred in the volume. A more flexible approach was 
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Figure 3. 

An example of a topological persistence measure defined on the superarcs of a contour tree or Reeb graph. The 
measure of persistence defined here uses the isovalue range—this is given as the difference of isovalues of the 
two critical vertices defining the superarc (see Table 5). 


used by Weber et al. [45] applying distinct transfer functions to each object, or 
topological zone, as directed by the user. This customization, implemented via the 
user interface, enables grouping of similar features and related components using 
colour and transparency. 


2.8 Temporal univariate scalar data 


Recent interest in topological visualisation research has focused upon the 
comparison of scalar data through topological methods to assign a degree of 
similarity to data. This has potential uses in different scientific domains including 
physics, chemistry, and climate science, which often feature a temporal component. 
The merge tree, with its simple data-structure, presents an attractive method for 
computing topological differences between data sets. Beketayev et al. [46] define 
a distance measure between merge trees with potential applications in a range 
of scientific disciplines. The algorithm uses a branch decomposition approach to 
deconstruct the merge tree into multiple sub-graphs, each a descending path from 
a saddle to a leaf vertex. Each branch decomposition can then be scored via an 
adapted form of the edit distance [47] between two graphs. A recursive algorithm 
then tests if two merge trees are within an epsilon value to determine if they are 
classified as similar. 

This approach was further extended by Saikia et al. [48] to produce the extended 
branch decomposition graph, a union of all individual sub-trees. This data structure 
allows quicker comparison between merge trees by computing multiple similar- 
ity thresholds in parallel. In addition the extended branch decomposition graph 
improves upon memory usage by reducing the redundancy found in multiple 
disjoint branch decompositions. Branch decomposition methods have also been 
applied to the more complex contour tree to compute similarity and symmetry in 
scalar topology [49]. 

The method used by Beketayev et al. [46] had a runtime of o(n”), where a merge 
tree contains n nodes, this was lowered to 0((n logn)’) in related work [48]. Potential 
downfalls to the optimised algorithm are that instabilities can arise from permuted 
forms of branch decompositions, this is handled in [46] by considering all possible 
permutations at the cost of scalability. The optimised algorithm uses the extended 
branch decomposition graph, a union of all possible branch decomposition graphs, to 
store all branch decompositions in a single tree structure. Two potential uses for the 
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Superarc Top supernode Bottom supernode 
Id Persistence Id Isovalue Id Isovalue 
1 0.67 12 0.84 13 0.17 

2 1.07 4 1.67 5 0.60 

3 0.54 6 1.25 7 0.71 

4 0.03 16 0.74 7 0.71 

5 0.11 7 0.71 5 0.60 

6 0.43 5 0.60 13 0.17 

Table 5. 


Persistence measures associated with Figure 3. 


distance measure are suggested; identification of similar structures within a single 
data set or, for time variate data, repetition between time steps. 

Fluid dynamics is an area of science that can benefit from visualisation of 
3D volumes with a time component. The complexity of mixing two fluids is one 
specific problem that can be better understood using topological methods as the 
system evolves over time [50]. Topological analysis makes it possible to take a 
slice of the volume at a given time-step and count the number of bubbles present. 
Alternatively, topological analysis enables tracing of properties, such the volume 
and centre of gravity of individual bubbles, throughout the simulation. 


3. Multivariate topological visualisation 


Multivariate topological visualisation is a more recent research interest in the 
visualisation community in comparison to topological visualisations of a single 
field. Due to the relative infancy of the field, Carr et al. [51] highlight the need for 
sufficiently complex data sets to extensively test the emerging algorithms in the 
area. 

Many of the topological structures and algorithms applicable to a single variable 
can be generalised to more than one variable being defined at each sampling point. 
We begin this section by defining the concept of Jacobi nodes, we then continue 
to examine structures used to capture and the display the topology of multiple 
variables. 


3.1 Jacobi node 


These are critical points in the scalar fields where the gradients of multiple 
inputs become parallel or equal to zero. A pre-requisite for the algorithms are that 
the fields are defined on a common sampling mesh meaning that the critical points 
are guaranteed to have the same geometric locations. 


3.2 The Reeb space 


The Reeb space is a generalisation of the Reeb graph to allow for multivariate or 
temporal data. The first discussion of using the Reeb space to compute topological 
structure of multiple functions is presented by Edelsbrunner et al. [52], where it is 
suggested that the Reeb space can be modelled mathematically in the formf:m + R*, 
where M represents the domain and f the output of k scalar functions. For the simple 
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case, where k = 1, this is directly comparable to the Reeb graph. The Reeb space 
extends this formulation to situations wherek > 2. 


3.3 The joint contour net 


Carr et al. [53] presented the first discrete representation of the Reeb space using 
the Joint Contour Net (JCN). For functions of n variables defined in an R” dimen- 
sional space the algorithm approximates the Reeb space as a number of multivariate 
contours named joint contour slabs. These represent connected regions of the 
domain with respect to the isovalue of multiple functions. In situations wheren = m 
the JCN can still be computed; however, the output is not an approximation of the 
Reeb space but instead a subdivision of the input geometry over n variables. The JCN 
captures the Reeb space as an undirected graph structure, where vertices represent 
slabs of n isovalue tuples, and edges are used to show adjacency between regions. An 
example JCN of two scalar functions is presented in Figures 4 and 5. 

In comparison to contour tree algorithms, the JCN is better suited to parallelisa- 
tion as each joint contour slab is constructed from smaller independent regions 
called fragments [54]. Existing development of the algorithm has focused largely 
on implementation as a parallel algorithm. A distributed memory model is used by 
Duke and Hosseini [55] to construct multiple sub-JCNs in parallel, these are merged 
into a single output as a final post-processing step. Current results suggest that the 
merge step is the limiting factor to parallel algorithm speed up; therefore, other 
parallelisation strategies are under investigation. 

In nuclear physics the JCN has previously been used to visualise and analyse 
scission datasets where it was used to identify the splitting of an atomic nucleus 
into multiple parts [56]. It was found that the JCN was well suited to capturing 
this divergent behaviour using proton and neutron density fields as inputs. This 
experiment was initially performed at a single temperature [57], but later repeated 


Figure 4. 

Two simple scalar functions are defined on a simplicial grid (left), where the dotted lines represent 
quantisation intervals. The quantised contour tree for each function (right) is shown mapped to the scalar field 
in the Centre. Image courtesy of Duke et al. [56]. 
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(6,8) (4,9) (3,9) (2,9) 


(4,1) : (5,2) (6,2) (7,2) (8,3) (9,3) 


Figure 5. 

A JCN capturing the bivariate topology of the two simple functions shown in Figure 4. The bivariate field is 
constructed by overlaying the quantisation intervals of the two input fields (dotted lines). A vertex is placed at 
the barycentre of each region, or joint contour slab, and edges mark adjacency. Image courtesy of Duke et al. 


[56]. 


at multiple temperatures [58] due to its ability to capture the splitting of the 
compound nucleus as a forking in the multi-field topology. Whilst performing the 
analysis a number of other events were captured and linked to the scission theory. 
More recently the JCN was used to visually analyse data from hurricane Isabel 
[59]. Vertices were used to represent the joint contour slabs by mapping to their 
barycentric spatial coordinates. An interactive environment was developed that 
allowed users to relate interactions in the temperature, pressure and precipitation 
fields to physical phenomena such as rain bands and the eye of the hurricane. The 
ability to relate properties of the JCN to known physical features helped to increase 
understanding of how the JCN is able to capture multi-field interactions. 


3.4 Related topological structures 


In multivariate topology a Jacobi set represents the set of critical points generated 
when one Morse function is restricted to the level set of another. Alternatively, the 
Jacobi set can be considered as the set of points where there is an alignment of the 
gradient of two or more Morse functions or the gradient of one function vanishes. 
The similarity between two or more functions can be evaluated using the resulting 
Jacobi set using a method defined by Edelsbrunner et al. [60]. Applications of the 
Jacobi set include use as a feedback loop on simulation parameters or to perform 
comparison of algorithms. An algorithm was presented by Edelsbrunner and Harer 
[61] for computing Jacobi sets of multiple Morse functions defined on a triangular 
mesh of isovalue m-tuples in Euclidean space. This allows multiple independent 
functions such as pressure, temperature, or wind speed to be used as inputs. 
Alternatively for temporal data multiple time-steps can be used as input to allow the 
evaluation of a function with respect to time. 

Carr et al. extend the JCN [53] to extract further topological structure from the 
Reeb space by first evaluating an intermediate structure, the Multi-Dimensional Reeb 
graph (MDRG) [62]. This is a hierarchical structure that recursively stores the joint 
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contours of each function (fofa ...,f,) restricted to the contours of those preceding it. 
At the top level the MDRG represents the contours of the function f} as a Reeb Graph; 
the second tier relates to the Reeb graph of function f, when restricted to the contours 
of f, continuing down to function f, restricted to the contours of fi, ...,f,1. Besides 
being used in the extraction of the Jacobi structure and the related Reeb skeleton [63], 
the MDRG represents a convenient structure for extracting the Reeb graphs of each 
individual function making up the Reeb space. 

An additional abstraction, the Jacobi structure, related to the mathematical 
topology concept of singular fibre-components is introduced by Chattopadhyay 
et al. [62] as part of the MDRG extraction algorithm. This represents an 
improvement over the Jacobi set by providing a method of relating the Jacobi 
set directly to the Reeb space. The Jacobi structure is able to capture the exact 
location of topological change and is defined as the projection of the Jacobi set 
from the domain to the Reeb space. In practice this extends the Jacobi set to also 
include the “regular sheets” connecting one another in the Reeb space. This 
means the Jacobi structure is able to capture elements of the topological struc- 
ture that the Jacobi set is unable to represent. The Jacobi structure is extracted as 
the set of critical nodes in the Multi-Dimensional Reeb Graph (MDRG), itself a 
structure for storing Reeb space criticalities. Forking in multi-field topology in 
nuclear scission data is an example behaviour that can be captured by the Jacobi 
structure [62]. 

The Layered Reeb graph is an alternative approach deployed by Strodthoff 
and Jüttler [64] for presenting the Reeb space of multiple scalar functions. This 
approach to representing the Reeb space differs from the MDRG [62] by work- 
ing directly with the Jacobi sets, rather than the more recently proposed Jacobi 
structure. 


3.5 Topological persistence and simplification 


Persistence in multivariate data sets is more complex to define in comparison to the 
univariate case. Simplification and persistence metrics can be defined on a number of 
secondary structures computable from the multi-field topology. The concept of isosur- 
face statistics [42, 43] is extended to multivariate inputs through the use of Continuous 
Scatterplots [65]. These can be defined to show relations between m dimensional inputs 
with n scalar fields; inthe case wherem = 3andn = 1the output approximates to the 
output of Meyer et al. [43]. 

Multivariate data gives rise to multi-filtrations due to their parametrisation 
by more than one variable; this leads to no definable compact invariants, such 
as the Betti numbers, existing in multi-fields. Therefore, existing concepts, 
such as the persistence bar code [66], do not directly generalise to multi-variate 
domains. However, this does not mean that persistence and simplification cannot 
be applied, instead other approaches have been suggested. The “rank invariant” 
is a method for representing persistence in a multi-field by generalising upon 
the concept of Betti numbers present in univariate topology. For the univariate 
case, the rank invariant and persistence bar code are the same [67]. The original 
algorithm used to compute the rank invariant was exponential in time complexity, 
this was later improved to polynomial time by reformulating the problem as an 
algebraic geometry problem [68]. 

The Jacobi set, where the gradient of multiple functions align or have a gradient 
of zero, can assist in defining persistence measures [60]. When the multi-field is 
used to represent temporal data this can be used to augment the univariate notion of 
persistence with a lifetime parameter. This approach was used by Bremer et al. [69] 
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to compute persistence in the context of the Morse-Smale complex. However, when 
generalised to non-temporal functions defining persistence as a feature of the Jacobi 
set becomes a non-trivial task [70]. 


3.6 The Reeb skeleton 


The Reeb skeleton (Figure 6) is a simplified graph structure that takes into 
account the size of connected components, allowing measures of persistence to be 
assigned to its arcs. An extended Jacobi set, the Jacobi Structure, is used by the Reeb 
skeleton algorithm to aid multivariate simplification [63]. 

The Jacobi structure [62] is a promising starting point for simplification due to 
its ability to separate the Reeb space, as approximated by the JCN, into singular and 
regular components. Just as in the univariate equivalent, the Reeb graph, singular 
nodes in the Jacobi structure map to topological changes in the multi-field. To 
exploit this, the Reeb skeleton extends the concept of the Jacobi structure further, 
primarily to aid multi-dimensional simplification [63]. 

The Reeb skeleton is generated as the dual graph of the singular and regular 
components captured in the Jacobi structure. Visually, this means the Reeb skel- 
eton translates the sheet-like form of the Jacobi structure into a simplified skeletal 
form. The simplified graph data structure allows measures of persistence to be 
assigned to arcs of the Reeb skeleton in a similar manner to that of the Reeb graph. 
Lip pruning techniques, similar to the leaf pruning method of simplification 
found in univariate topological structures [12] can then be applied to progres- 
sively remove noisy features in the multi-field. Example persistence measures that 
can be applied to the JCN include the accumulated volume of joint contour slabs 
in a connected region. 
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Figure 6. 

An example of a Reeb skeleton, showing how it relates to the JCN (a); (b) shows the full, non-simplified Reeb 
skeleton; (c) modifies the Reeb skeleton to highlight only critical changes in the multivariate topology (red 
vertices); (d) shows how performing simplification using the Reeb skeleton removes less significant regions of 
the topology; (e) after simplification, the Reeb skeleton can be reduced to only key vertices; (f) the arrowed 
regions relate to the two surfaces. Image courtesy of Chattopadhyay et al. [63]. 
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4. Visualisation in higher dimensions 


When moving to higher dimensional spaces it is beneficial to generalise the 
terminology used to describe the geometry. The n-dimensional analogue of the square 
is the hypercube, often shortened to n-cube (see Table 6). When working in higher 
dimensional spaces it is common to perform a simplicial sub-division of the n-cube 
into n-simplexes, this helps to avoid the ambiguities often associated with MC style 
algorithms [71]. 


4.1 Projection and perception 


Creating easily perceivable and topologically correct three-dimensional models 
for projection on to flat surfaces, the computer monitor, is a difficult task. Volumes 
of data become increasingly hard to visualise as their dimensionality increase; for 
example, by introducing a temporal fourth-dimension. One of the major limiting 
factors is our in ability to perceive things four-dimensionally. A metaphor for trying 
to understand four-dimensions in a three-dimensional world is to consider the case 
of two-dimensional creatures trying to understand a three-dimensional world. This 
is a thought exercise discussed in [72], which also discusses how a four-dimensional 
Euclidean representation of space-time relates to real world physics. 

Existing projection methods are available that take a four-dimensional objects and 
display them on a two-dimensional surface, usually in wire-frame form. Quite often 
the projections are animated to show the object as it rotates on one or more axis; how- 
ever, this can also be adapted to allow the user to rotate the object through conventional 
approaches such as mouse interaction. The effect of perspective and isometric projec- 
tion are explored by Hollasch [73] using tesseracts—the 4D hyper-cube. In addition, 
the use of ray-tracing in four-dimensions can produce understandable images, as depth 
cueing is handled automatically by the algorithm in the creation of shadows. However, 
the added complexity of a fourth dimension in a looping ray tracing algorithm makes it 
time consuming and potentially unworkable for real time display. 


4.2 Computing surfaces in higher dimensions 


Whilst presenting difficulties with regard to visualisation, it can be beneficial to 
compute surfaces in higher dimensions. This is especially true in the case of volu- 
metric data with a temporal element; animation can often be used to reconstruct this 
form of data but that presents its own perceptual issues. Computation of isosurfaces 
on fields existing in r* space can help to improve animations by providing interpola- 
tion between discrete time-steps, resulting in smoother and easier to perceive anima- 
tions. Alternatively the high-dimensional topology can be sliced along arbitrary axes 
to provide a snap-shot with a reduction in dimension (e.g. rR‘ + R°). 


Name Dimension n Common name Vertices Edges Faces Cells 

0-cube 0 Point 1 — — — 

1-cube 1 Line-segment 2 1 — — 

2-cube 2 Square 4 4 1 — 

3-cube 3 Cube 8 12 6 1 

4-cube 4 Tesseract 16 32 24 8 
Table 6. 


Properties of hypercubes in dimensions zero to four. 
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n Centroid division (2”'n!) n! 

2 4 2 

3 24 6 

4 192 24 

5 1920 120 

6 23,040 720 
Table 7. 


Number of simplexes generated using differing sub-division techniques in n dimensions. 


An unavoidable consequence of upping the dimensionality of the input field is 
an increase in the complexity of its storage and computation. An early example of 
computing surfaces in > 3 dimensions is considered by Weigle and Banks [74] using 
a recursive technique to split n-dimensional cells into n-simplexes. The resulting 
surfaces exist in four dimensions and are able to be displayed using two techniques; 
stereographic projection (Section 4.1) or slicing to reduce dimensionality. It is noted 
that the centroid division technique used to break cells into simplexes used in this 
work is suboptimal (see Table 7), but can be improved using lookup tables similar 
to those used in marching cubes. The technique was used in [75] to compute the 
swept volumes of time-varying data generated from electromagnetic field simula- 
tions, allowing them to be displayed as animations. 

Extending MC into the fourth dimension, and beyond, was explored by Bhaniramka 
et al. [76]. At the time of the report, a relatively small amount of work had been 
conducted in studying isosurfaces beyond the third dimension. It was found that by 
extending MC into 4 dimensions the look-up table, following removal of symmetrical 
configurations, required 222 separate configurations. As with the 3D variation of MC, 
one aspect to be taken into consideration is the dealing of ambiguous configurations; a 
mathematical proof of correctness is provided to verify that the topological structures 
generated are valid. An example use of this algorithm would be to present volume data 
with a time dimension by selecting three-dimensional slices of the hyper-volume. 


4.3 Topological landscapes 
An alternative approach of viewing contour information, the topological 


landscape, was proposed by Weber et al. [45], using the contour tree to build a 3D 
terrain model (Figure 7). The purpose of this is to harness the natural ability that 


Figure 7. 
Peaks within a topological landscape (left) correspond to distinct topological features in the direct volume 
render (right). Image courtesy of Weber et al. [45]. 
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humans have in understanding terrain structure and use it to provide an easier to 
understand model of the underlying data topology. Valleys in the terrain illustrate 
events where a contour splits into two or more parts and peaks represent where two 
or more contours merge. Topological landscapes can be applied to contour trees of 
any number of dimensions; hence, they can provide a useful method for exploring 
what is happening in high dimensional data sets. The topological landscape meth- 
odology was further expanded using a number of different methods for laying out 
the data, primarily from established 2D visualisation techniques [77]. 


5. Conclusion 


In this chapter we first considered problems existing in marching cubes 
algorithms relating to the topological correctness of the data. We demonstrated 
how, through the use of topology, a correct representation of a scalar field can be 
captured using graph structures. These could be used to seed topologically correct 
meshes for rendering, or to provide a means to analyse the data. After providing a 
description of algorithms for data sets consisting of a single variable, we showed 
how many of the techniques can be generalised to multivariate data. Finally, we 
considered the techniques for displaying data existing in higher dimensions. 
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